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Abstract 

Predictions of a thermodynamical model of hadron production for multiplicity 
distributions in e + e~ annihilations at LEP and PEP-PETRA centre of mass ener- 
gies are shown. The production process is described as a two-step process in which 
primary hadrons emitted from the thermal source decay into final observable parti- 
cles. The final charged tracks multiplicity distributions turn out to be of Negative 
Binomial type and are in quite good agreement with experimental observations. 
The average number of clans calculated from fitted Negative Binomial coincides 
with the average number of primary hadrons predicted by the thermodynamical 
model, suggesting that clans should be identified with primary hadrons. 
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1 Introduction 



Multiplicity Distributions (MDs) are fundamental observables in multiparticle production in 
e + e~ collisions. Theoretical investigations of their properties, in the QCD framework, are based 
on a comparison of approximate calculations at parton level with experimentally observed fi- 
nal particles distributions via the assumption of the Local Parton Hadron Duality [|TJ or its 
generalization 0. In this context the important role of resonances and particles decay is not 
explicitly taken into account. This fact is particularly unsatisfactory because a large fraction 
of final particles are indeed decay products of heavier ones. A second approach is based on 
Monte Carlo models, like JETSET [Q and HERWIG Q, which account for the decay chain 
following the perturbative phase and hadronization. However, predictions of these models are 
obtained by using hadronization schemes requiring a large number of free parameters in order 
to reproduce experimental data. 

In this paper we discuss the MDs obtained in a thermodynamical model of multiparticle produc- 
tion in e + e~ annihilations, introduced by one of the authors ||, which successfully reproduces 
the production rates of the various hadrons species both at PEP-PETRA and LEP centre of 
mass energies. The model is based on the identification of jets with thermalized hadron gas 
phases after the hadronization of primary quarks has taken place. As a first approximation it 
is assumed that only two phases with opposite momenta are generated in an e + e _ — > qq event, 
namely multi-jets events are neglected. The observed multi-particle production in hadronic 
events is the result of a two-step process: in the first one some primary hadrons (particles 
and resonances) emerge directly from the thermalized phases after having decoupled. These 
primary hadrons then decay according to known decay modes and branching ratios, giving rise 
to observable particles in the detector. Only three parameters describe those phases at the 
decoupling time: the temperature T, the volume V and a strangeness chemical equilibrium 
suppression parameter j s . They are determined by fitting the calculated average production 
rates of each hadron to those measured both at LEP (y^i ~ 91 GeV) and at PEP-PETRA 
(-^/s = 29 -T- 35 GeV). A major role in the determination of the hadron rates is played by the 
conservation laws: the quantum numbers of the jet are related to those of the primary quark 
from which the jet itself originated. Accordingly, it is assumed that each jet keeps the charm 
and beauty of the parent quark while non- vanishing baryon number and strangeness are allowed 
provided that the baryon number and strangeness of the whole system are zero. 

2 Thermal fluctuations and multiplicity distributions 

Consider a hadron gas at temperature T and volume V. If T « (9(100) MeV, i.e. less than 
the mass of all hadrons but pions, the simple Boltzmann statistics holds for all hadrons except 
pions, which indeed obey Bose statistics ||. Therefore, if the quantum numbers of the gas are 
not fixed, each species of hadron fluctuates independently according to a Poisson distribution, 
whereas pions fluctuate according to a different distribution that we define as the 7r-distribution 
f n and we derive in Appendix A. On the other hand, if the flavour quantum numbers of the 
gas are fixed, conservation laws generate correlations between different hadron species, so that 
fluctuations are no longer independent. 

Therefore, the probability to observe n\ tt + , n 2 tt°, n 3 tt~, n 4 hadrons of kind 4,..., uk hadrons 
of kind K in a hadron gas system is given by: 
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P{m,...,n K ) = ^QyAK^2,n 3 )n ^r e Zi 5 Q,j2l in ^ > W 

where 



* = (2Jj + 1) j— - / e — , (2) 



(2tt) 3 

Q = (N, S, C, B) is a four dimensional vector with integer components representing baryon 
number, strangeness, charm and beauty of the gas; q^, Jj, are the quantum numbers vector, 
the spin and the mass of the i th particle and f n is the 7r-distribution. The factor 5^ ^ K n q . ac_ 

counts for the conservation of the mentioned quantum numbers. Z(Q) is the partition function 
of the system, namely: 



k m 



z(Q) = E ••• E AKn 2 ,n 3 )n fy^^E .««* • ( 3 ) 

With the substitution 



ni=0 n x =0 i=4 



5 QE>^ = (2^)4 / d ^ ex P^(Q - E^^) • 4>)} , ( 4 ) 
where = (0i, 02, 03, 04), Eq. (3) becomes: 

Z(Q) = (2^/ ^ ^exp{E^e"^ n , (5) 

which is the expression obtained for the hadron gas partition function in [BJ starting from the 
general formulae of partition functions of thermodynamical systems with internal symmetry ||. 
The function turns out to be: 



3 V f \/p*+mf 

F n = exp{-E^p J d 3 p log (1 - e — )} ; (6) 

where the sum runs over the three pion states. 

From Eq. (1) one can build up the expression of the joint probability to observe a i^-uple 
n = (ni, . . . , rix) of numbers of primary hadrons in the first jet and a .fT-uple m = (mi, . . . , mx) 
in the second jet 

K n< K rm 

P(n,m) = ^/.(r2i,n 2 ,n3)^(m 1 ,m 2 ,m 3 )n -^e"** ]J -^e~ Zi 5 0)S n . q . +m . q . , (7) 
with Z given by 
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E ••• E E ••• E U{ n ^ n 2,n 3 )U{m 1 ,m2,'m 3 ) 

ni=0 iij(=0mi=0 raj(=0 



if n, if m 



i=4 n *' i=4 m ' - 



It should be pointed out that the partition function of the two-jet system is further modified 
by other constraints: 



- net baryon number per jet | N | < 2 

- net strangeness per jet | S \ < 3 

- S 7^ in e + e~ — > ss events 

- net charm per jet | C \ = 1 in e + e~ — > cc, \ C \ = otherwise 

- net beauty per jet | B \ — 1 in e + e~ — > bb events, | B \ — otherwise 

These constraints can be easily implemented by multiplying the right-hand side of Eq. (7) by 
additional suitable 5's and by recalculating accordingly the Z in the Eq. (8). Furthermore, the 
probability of production of primary hadrons containing n strange valence quarks should be 
multiplied by a suppression factor 7™. 

The Eq. (7) suggests a simple procedure for calculating the charged tracks MDs in e + e~ — > qq 
events with a Monte-Carlo method. The first step of the procedure is to pick up two if-uples 
((ni, . . . , n K ), (mi, . . . , m K )) randomly according to Poisson or 7r-distribution and to accept the 
event only if the conditions imposed by the 5's in Eq. (7) and by the other requirements are 
fulfilled. The second step is to perform all decays of the generated primary hadrons according 
to known branching ratios, until 7T, K, K®, fi or stable particles are reached, in order to match 
the MDs measured by e + e~ colliders experiments, including all decay products of particles with 
ct < 10 cm. 



3 Results 

We calculated the multiplicity distributions by the Monte-Carlo method described in the previ- 
ous section at LEP and PEP-PETRA centre of mass energies. We used the same values of the 
parameters T,V and j s determined in |7] (see Table 1) by fitting the calculated average produc- 
tion rates of various hadron species to the measured ones ||. Two additional parameters are 
introduced in order to take into account the possibility of different values of V as a function of 
the primary quark mass, i.e. two dimensionless variables x c , Xb G [0, 1] such that x c V is the vol- 
ume in e + e~ — > cc and x^V in e + e~ — > bb events, whereas V(l— x c R c — x^R^) / (R u +Rd+R s ) is the 
volume for e + e~ — > qq where q is a light quark and R q = a(e + e~ — > qq)/a(e + e~ — > hadrons). 
The factor x c (xb) has been determined constraining the difference 5 cq (5b q ) between the av- 
erage charged multiplicity in e + e~ — > cc (e + e~ — > bb) and the average charged multiplicity 
in e + e~ — > qq, with q = u,d,s for LEP, q — u,d, s, c, b for PEP-PETRA to be equal to the 
measured ones P,p!0|. The measured values of 5 cq and Sb q have been averaged according to the 



procedure described in [ |TT|j . 

We generated 10000 events for each centre of mass energy in order to get the MDs of primary 
hadrons and of charged tracks both in single hemisphere (i.e. single jet according to the ther- 
modynamical model) and in full phase space. 

The calculated MDs are fitted with a Negative Binomial (NB) distribution with Maximum 
Likelihood method (see Fig. la, lb and 2a, 2b). Average charged tracks multiplicity n and 
parameter k (linked to dispersion D by the relation D 2 = n + n 2 /k) obtained from the fit are 
quoted in Table 2 together with corresponding values from analogous NB fits to experimental 
data j |T2| , |H| | . Two kinds of error affect the determination of the NB parameters: the first one 
is the fit error, due to the limited statistics of 10000 events; the second one is related to the 
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uncertainty on the parameters of the thermodynamical model T, V, 7 S , x C} Xf, used as input in 
the event generation. This latter error, a systematic one, has been estimated by varying the 
values of parameters by their error, as quoted in Table 1, generating a new sample of events 
and repeating the NB fit. The covariance matrix M of n and k has been determined with the 
usual formula 

M = JCJ T (9) 

where C is the covariance matrix of T, V, j 3 , x c , xj, and J is the jacobian matrix relating (n, k) 
to (T, V, 7 S , x c , x b ). 

The obtained systematic errors dominate over statistical ones, as shown in Table 2. Since the 
X 2 test on the NB fit consistency is also affected by those systematic effects, we estimated 
the uncertainty on x 2 according to the same procedure used for n and k. Taking this error 
into account, the x 2 test indicates a clear compatibility with the NB distribution within the 
uncertainty of the input parameters. 

The agreement between the parameters obtained by the model and those obtained from NB 
fits to experimental data is remarkable (see Table 2). However, discrepancies arise in the direct 
comparison between calculated distributions and measured ones as shown in Fig. 3. The 

overestimation of the probability distribution in the low-multiplicity bins at PETRA energies 
and the consequent increase of the dispersion D, can be attributed to violations of the statistical 
framework of the canonical ensemble at low centre of mass energies. It is expected indeed that 
the probability of producing a low number of primary hadrons (for instance or 1) should be 
strongly suppressed in comparison to the predictions of the canonical distribution, due to the 
exact energy conservation. This deviation from the canonical behaviour is much less visible at 
LEP energies, where the average number of primary hadrons is higher. 

Another discrepancy is related to the so-called "shoulder effect": the experimental charged 
tracks MD deviates from the NB distribution due to a shoulder structure in the high multiplicity 
region, which is interpreted as the effect of the superposition of low-multiplicities 2-jets events 
with multi-jets events yielding higher multiplicities. On the other hand, it has been shown 
that good NB fits are obtained with MDs in selected 2, 3 or 4 jets event samples as well as 
with individual jets MDs [fLql . Since the thermodynamical model is based on a 2-jet scheme 
and neglects multi-jets events, an analogous disagreement with data is expected, while the 
agreement with NB confirms its validity in reproducing MDs in 2-jets events. It should be 
mentioned also that a comparison between our predictions and a 2-jet events data sample is 
not possible because the parameters of the model have been tuned to reproduce the production 
rates of the hadron species in the overall e + e~ — > qq data sample. 



4 Comparison with the clan model 

The clan model has been introduced in order to interpret the wide occurrence of the NB 
regularity of charged tracks MDs in all high energy reactions jn|. This model is based on 
the assumption that initial independent, i.e. Poissonian, production of primary objects called 
clans, is followed in the second step by their decay into final particles according to a logarithmic 
distribution. The average number of clans is determined in terms of the n and k parameters of 
the observed Negative Binomial distribution according to the formula: 

N c = klog (1 + n/k) . (10) 

Since the thermodynamical model also contains a two step structure of the hadron production 
process consisting in an initial generation of a number of primary hadrons from the thermal 
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source which then decay into final observable particles, one can ask whether clans can be 
identified with the primary hadrons of the thermodynamical model. The average number of 
primaries, in this case, is determined independently by using the parameters of the thermody- 
namical model determined as mentioned in Sect. 3. 

Results of this comparative analisys are summarized in Table 3. The agreement between the 
average number of primaries as resulting from the two models is striking. It should be noticed 
that in the clan model primaries are distributed according to an exact Poisson distribution, 
while in the thermodynamical model the distribution of primaries can be fitted by a NB dis- 
tribution with a quite large k value, which can be well approximated by a Poisson distribution 
(see Fig. lc, 2c). This small deviation from an exact Poissonian behaviour is a consequence of 
the presence of conservation laws as specified in the 5's of Eq. 8, and of correlations contained 
in the ^-distribution for primary pions. It should be noticed, however, that the ^-distribution, 
with the actual values of T and V, is very close to a Poisson distribution (see Table 4). 

5 Conclusions 

The thermodynamical model of hadron production in e + e~ collisions has been applied to the 
study of the charged tracks multiplicity distributions. Predictions of the model reproduce the 
gross features of experimental data both at LEP and PEP-PETRA centre of mass energies. The 
observed discrepancy with experimental measured MD at LEP energy is due to the presence of 
multi-jets events which are neglected in the model (shoulder effect). A natural violation of the 
canonical distribution, expected at low energy, is responsible for the observed deviation from 
the experimental measured MD at PETRA in the low multiplicity region. Also, the predicted 
charged tracks MDs are shown to be in very good agreement with the long-standing observed 
NB regularity and both at LEP and PEP-PETRA centre of mass energies. The primary hadrons 
MD is also well fitted by a NB distribution with a quite large value of k parameter, close to the 
Poissonian limit. The average number of primary hadrons, independently determined in the 
thermodynamical model, is approximately equal to the average number of clans calculated from 
the parameters of the NB fit to the final charged tracks MDs. This identity suggests identifying 
clans with primary hadrons and confirms the intuition that the occurrence of Negative Binomial 
distribution is deeply related to a two-step process ||T6| . 

These results confirm previous successful predictions of the thermodynamical model, opening a 
new perspective in the understanding of the hadronization process at high energy. An extension 
of this model to include multi-jets events is going to be pursued in order to reproduce also fine 
features of experimental multiplicity distributions. 
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Appendix 



A Determination of the ^-distribution 

In a boson gas at temperature T and volume V the probability to observe rt\ particles in the 
kinematical state 1 of energy e 1 , n 2 particles in the kinematical state 2 of energy e 2 , ... is: 

P( ni , . . . , n M ) = H(l - e-^' T ) e- n ^' T . (11) 

i 

The probability that the overall population is N can be written as: 

oo oo 

P(N) = E ■ ■ ■ E IK 1 - e~ £ > /T ) e- n ^ /T S N , E ^ • (12) 

ni=0 n M =0 i 

Since: 

1 f 27T 

VE t «, = 2^J d< ^ expi-KN - E n i)^i > ( 13 ) 

Eq. (12) becomes: 

P(N) = ± / d<j> e- iN<t > TT ; • (14) 

In the limit of a continuum of energy levels Eq. (14) can be rewritten as: 

p(N ^hr d ^" N * e ^ 2j+i ^i d3p iog i i ^^ } - <i5) 

where e = \Jp 2 + m 2 and J is the spin of the particle. 

Eq. (15) in this form can be transformed in an integral over the unitary circle in the complex 
plane. Let w = exp{i0} and 

J = exp{(2J + l)^p Jd 3 p log(l-e- £ / T )}. (16) 

Then: 

P(iV) = J^ 6XP{ - (2J+ 1} (2^ / ^ bg(1 ~ e ' £/T ^ ■ ( 17 ) 

The residuals theorem can be used to get the final expression of P(N): 

I d N V r 

It turns out that the derivatives in w = can be expressed as a function of the functions z^ n ) 
defined as: 

V r „./^T=Z„ _ VT 



Z{n) (2tt)3 



/ d* P e-V^/r = ( 2 J + 1) XL m*K 2 C^) (19) 



where K 2 is the modified Bessel function of order 2. 

Finally, by using the above expression of P(N) one gets the 7r-distribution fn(ni,n 2 ,ri3) 



P+(ni) Po(n 2 ) P_(n 3 ). In Table 4 we show probabilities P+(N) of producing N ir + up to N = 6 
at T = 163 MeV and V = 20 Fm 3 . The probability values are compared with those obtained 
from a Poisson distribution with the same average value N. 
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Figure captions 



Figure 1 Multiplicity distributions calculated in the framework of the thermodynamical model at 
y/s = 91.2 GeV (dots) and corresponding Negative Binomial fits (solid lines): a) charged 
tracks in full phase space, b) charged tracks in a single hemisphere, c) primary hadrons. 

Figure 2 Multiplicity distributions calculated in the framework of the thermodynamical model at 
■s/s = 29 -T- 35 GeV (dots) and corresponding Negative Binomial fits (solid lines): a) 
charged tracks in full phase space, b) charged tracks in a single hemisphere, c) primary 
hadrons. 



Figure 3 Charged tracks multiplicity distributions in full phase space; comparison between predic- 
tions of the thermodynamical model (histogram) and data (dots) [[13,14] at a) ^/s = 91.2 
GeV and b) = 29 + 35 GeV. 
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Tables 



Parameters 


v^ = 91.2 GeV 


y/s = 29 4- 35 GeV 


Temperature(MeV) 


162.9±2.1 


169.3±3.5 


Volume(Fm 3 ) 


21.4±1.9 


9.3±1.4 


7* 


0.696±0.027 


0.811±0.046 


x c 


0.884±0.029 


0.89±0.12 


Xb 


0.695±0.016 


0.47±0.15 



Table 1: Values of the parameters of the thermodynamical model [[?]]. 



Full phase space 




y/s = 91.2 GeV 


y/s = 34 GeV 


Parameters 


Calculated 


Measured 


Calculated 


Measured 


n 


21.45 ±0.093 ±0.58 


20.71 ±0.77 


12.79 ±0.08 ±0.50 


13.59 ±0.46 


k 


28.34 ±1.43 ±1.41 


24.33 ±0.71 


13.03 ±0.58 ±0.66 


52.63 ±5.6 


D 


6.14 ±0.07 ±0.17 


6.28 ±0.43 


5.04 ±0.06 ±0.50 


4.14 ±0.39 


X 2 /dof 


1.81 ±0.60 


2.35 


10.9 ±5.53 


2.38 


Single hemisphere 




y/s = 91.2 GeV 


yfs = 34 GeV 


Parameters 


Calculated 


Measured 


Calculated 


Measured 


n 


10.72 ±0.045 ±0.29 


10.35 ±0.47 


6.39 ±0.035 ±0.25 


6.78 ±0.33 


k 


16.42 ±0.65 ±0.45 


15.06 ±0.39 


9.14 ±0.35 ±0.78 


50.0 ±5.0 


D 


4.21 ±0.035 ±0.056 


4.19 ±0.32 


3.30 ±0.03 ±0.23 


2.78 ±0.22 


X 2 /dof 


1.76 ±0.30 


2.87 


8.83 ±2.8 


4.68 



Table 2: Results of Negative Binomial fits to calculated charged tracks multiplicity distributions 
compared to results of analogous fits to experimental data [12|, [13|]. The first quoted error is 
the fit error, the second one is due to the uncertainty on parameters of the thermodynamical 
model. The dispersion D is determined with the formula D 2 = n + n 2 /k. 
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= 91.2 GeV 


y/s = 29 4- 35 GeV 




15.97 ±0.20 ±0.49 


8.91 ±0.12 ±0.78 


N P 


17.10 ±0.07 ±0.55 


9.39 ±0.05 ±0.36 


k 


41.42 ±3.09 ±9.45 


32.18 ±3.11 ±26.0 


X 2 /dof 


1.63 ±0.60 


2.04 ±0.99 



Table 3: Average number of clans N c extracted from the Negative Binomial fits to charged 
tracks multiplicity distributions compared with the average number of primary hadrons N p . 
Also shown are the results of a Negative Binomial fit to the calculated primary hadron distribu- 
tion (k and x 2 )- The first quoted error is the fit error, the second one is due to the uncertainty 
on parameters of the thermodynamical model. 



N 


P + (N) 


Poisson 





0.311 


0.293 


1 


0.346 


0.359 


2 


0.208 


0.221 


3 


0.0903 


0.0904 


4 


0.0320 


0.0278 


5 


0.00998 


0.00682 


6 


0.00287 


0.00139 



Table 4: Probability of producing N ir + at T = 163 MeV and V = 20 Fm 3 compared to a 
Poisson distribution with the same average value. 
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Figure 1: 
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Figure 2: 



14 




Figure 3: 
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